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The continuously changing ionosphere profile is approximated by a series of homogeneous layers, 
and the radial field variation in the individual layers is obtained by solving a fourth order equation. 
Its biquadratic solutions are adequate only in regions where TE and TM modes are coupled, but more 
accurate solutions are required in lower layers of the ionosphere. The field components below the 
ionosphere are computed after multiplying a sequence of 4 X 4 matrices, each of which represents the 
effects of one of the ionosphere layers. However, simpler multiplications can be used for heights 
above the D layer. After computing the impedance matrix at the lower boundary of the ionosphere, 
the propagation parameters are computed with neglected coupling between TE and TM waves below 
the ionosphere. 

The radial magnetic field has the most pronounced effects at nighttime in the lower end of ELF 
range, where it increases the attenuation significantly, and also in the VLF range, where it decreases 
the attenuation. At daytime the presence of the radial magnetic field has only minor effects on wave 
propagation. 

1. Introduction 

The ionospheric anisotropy has been investigated for a model of a homogeneous ionosphere 
by Schumann [1955] and Wait [1962, 1963a]. The radial magnetic field causes a coupling between 
the TE and TM modes of propagation, but this coupling has a small effect for frequencies in the 
ELF and VLF range. Volland [1962, 1964] considers propagation below an inhomogeneous and 
anisotropic ionosphere using matrices of reflection coefficients. A formalism for treating a 
multilayer ionosphere has been presented by Wait [1963b], and the formal solution for an M 
layer ionosphere is shown to require the solution of 4M + 2 simultaneous equations in 4A/ + 2 
unknowns. Another solution intended for frequencies in the Schumann resonance range (30 c/s 
and lower) has been derived in the thesis of Thompson [1963 J, and shown by Galejs [1964a]. His 
procedure considers each resonant mode separately and requires the multiplication of a sequence 
of 4 X 4 matrices. This solution does not consider the displacement currents of the lower iono- 
sphere regions and results in (^-factors of the Schumann resonances that are larger than observed. 

In the present paper the effects of a radial magnetic field are considered with a multilayer iono- 
sphere model. Matrix multiplication techniques are used to relate the field components at the 
lower ionosphere boundary to the wave amplitudes of the upper ionosphere layer. The propa- 
gation parameters of the waves in the VLF and ELF frequency bands are obtained after computing 
the surface impedance at the ionospheric boundary (allowing for displacement currents) and by 
using the appropriate modal equation in an iteration procedure. It is shown that the biquadratic 
solution of two coupled equations or of a single quartic equation [Wait, 1963b; Thompson, 1963; 
Galejs, 1964a] is not sufficient for determining the propagation parameters of the elementary 
waves in the lower ionosphere regions, where the TE and TM waves tend to become uncoupled 
by the local medium. These solutions are described in section 2. In applying the standard 
technique of matrix multiplication and inversion to the upper layers of the ionosphere model, 
excessive cumulative errors occur particularly for frequencies in the VLF range. However, 
for these parameters the matrix multiplication and inversion can be replaced by a simple one 
term multiplication of the matrix elements, which in addition to speeding up the computational 
processes avoids the above convergence problem, as indicated in section 3. Numerical examples 
of propagation parameters for the ELF and VLF bands are shown in section 4. 
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2. Fields in the Anisotropic Ionosphere 

The wave solutions for the fields in the anisotropic ionosphere can be derived from Maxwell's 



equations: 



and 



VxE = icofjioH 
VxU=-ioje [e]E. 



(1) 

(2) 



For an assumed exp(— icot) time dependence of the fields and an r-directed static magnetic 
field the permittivity tensor [e] is given by the classical magneto-ionic theory as 
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The conductivity components cr , <Xi, and a 2 of a typical lower ionosphere model have been reported 
by Galejs and Row [1964] and Galejs [1964b]. The field components E r and H r are proportional 
to fr{6) = Pv(— cos 0), while E$, He, E^, and H^ are proportional to f$(0) =f<f>(d) = dP v {— cos)/d0. 
After writing out the 2 vector equations (1) and (2) in their 6 scalar components, the 6 functions can 
be canceled leaving 6 scalar equations in the r-dependent components EJ(r) = Ejlf/fl), HJ(r) — Hjlfj(d), 
where j = r, 6, or (/>. This gives 
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After introducing the notation 
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0=1, 2, 3), a manipulation of (4) to (9) results in 
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2?$ can be eliminated from (11) to (13), which gives a fourth order differential equation for computing 
E°. The parameter v is defined by v(y+ l) = (k aS) 2 , where k = a> V/x<)€o, a is the radius of earth, 
and S is the yet undetermined root of the modal equation. The factors v(v-\-Y)j{kjr) 2 change 
gradually relative to {r 2 E () t ), and only the first order change is considered. These factors and their 
first derivatives are replaced by constants 



Qj = 



v{v+l) 

' (kjrm) 2 ' 



(14) 



Qj 



Tr^'^ 



(15) 



where r m is the mean radius of an ionospheric layer. This results in 

^ (r*£») + b, £ (r 2 £») + B 2 ^ (r 2 £',!) + b, j- (r 2 £») + (B„ + KV'E" = 
dr dr 1 dr z dr 
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A first order solution of (16) is obtained by neglecting the small terms proportional to b fi . This 
solution is of the form 
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Equation (23) is seen to be in agreement with (34) to (36) of Wait [1963]. For negligible displace- 
ment currents kj » k = aj Vjito€o. Noting that v(v-\-\) ~ (k r) 2 , qj << 1 and (23) simplifies fur- 
ther to the form given by Thompson [1963] . 

The second order solution (designated by primes) is obtained by a single iteration substituting 
(22) in (16), solving for k 2 , and by using (23) for k in the right-hand side of this equation. This 
procedure gives 



k; 2 =^(k]-ib,k^ib i k j ^B {) ^bo) = k 2 ^(b ] -b,kj)^^ 



(24) 



707 



where 7=+ or—. The second solution which corresponds to the negative of kj is designated as k" 2 . 
Although there are only minor differences between the solutions (23) and (24), it is essential to use 
(24) in the subsequent expressions (35) and (48), because (23) does not result in decoupled TE and 
TM field components as e 2 — > and ei and 63 — * 1 in (3). 
In the limiting case of fe— »0, (12) simplifies to 



-^(r 2 E°) + kj 



^IvrV^ (25) 



and, after applying (4), (11) gives 
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In the planar approximation (25) and (26) have solutions 

r 2 E° r = e ikr (27) 

with 

k* = k%l-q 3 ) (28) 

and 

r*H%=e*r (29) 

with 

k 2 =k 2 a- qi ). (30) 

The same solutions, (28) and (30), are also obtained for fe - *0 from (23). The upper and lower 
signs of (23) apply therefore to TM and TE field components derived from E* and //£, respectively, 
which become coupled for fe ^ 0. 

The field components will be derived using the planar approximation, (14) and (15). The field 
components will be related to the amplitudes of the radial electric and magnetic field components. 

The component E r is assumed to be of the form 

E« = a v v(y+\)-e* kr , (31) 

where k has four possible values according to (23). H% follows from (7) as 

HI = a u i(oeoe 3 - e ikr , (32) 



E® d is obtained from (6) 
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Computing the second derivative of (r//g) from (9), applying (4) and substituting (5) and (6) for the 
derivatives of (rE$) and (rE§) gives 

F 

H° d = a v io)e € 3 — e lkr (34) 
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with 



&K1- g3 ) 
[JP-AJ(l-flO]' (35) 

£$ and //J! follow from (5) and (4), respectively as 
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r o)/jL kr z y ' 

When using (24) for k 2 in (35), there is no coupling between the TM components (E% Z?jj, //$) and the 
TE components (E%, H% ffg) for e 2 or fcf-» in (4) to (9), and £J f #? and H° d of (36), (37), and (34) 
approach zero. 1 For negligible displacement currents, qj<< 1 and the expressions (31) to (37) 
are simplified to 

E° = aMv+l)-r 2 e ikr (38) 

Eg = a„-^^e* r (39) 

E% = a v -^7-^e*' (40) 

e 2 ±iei r 
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r 

These expressions of circularly polarized field components have been obtained also by Thompson 
[1963]. 

In an alternate representation, H r is assumed to be of the form 

H» = b v v{v+ l)jge ikr ; (44) 

E% and H° follow from (4) and (5) as 



E%=-b/-^e*r (45) 

W = b v je*>: (46) 



1 For decoupled field components, F of (35) should approach zero. However, the field components are not decoupled when using (23) for k 2 in (35). For k :i = k t , 
F is finite and independent of k 2 . For k :i ^ k } and kt~* 0, F ~ k%l(k% — k 2 ). After relating kj to the tensor conductivity components (Tj by (10) and (3), kj can be ex- 
pressed in terms of collision frequencies v, plasma frequencies to and gyrofrequencies um of the charged particles of the ionosphere, as shown by (3) to (7) of Galejs 
and Row [1964]. The last expression given above for F becomes F — v/\(dh\ which is increasing with a decreasing altitude contrary to expectations. 
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Computing the second derivative of (rH%) from (8) and applying (6), (7), and (38) gives 

b v v(v + 1) cofipk E 



where 



H% and E° d follow from (7) and (6) as 
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When using (24) for k 2 in (48), there is no coupling between the TE and TM field components if 
€2 or A|— > 0, and 2?°, EJJ, and //£ approach zero. For negligible displacement currents, q$ < 1 and 
the expressions (47) to (50) can be simplified to 

Eo =± bMv+lUo>f«Jc eikr (5i) 

£|=+^e* (52) 

fl>=^e*. (53) 

By noting that 

6„ = ± a„A|/ ( ia>/xo&) (54) 

the field components (38) to (43) are seen to be identical to (44) to (46) and (51) to (54). 



3. Surface Impedance of the Stratified Ionosphere 

The ionosphere is subdivided in a number of concentric spherical shells as is indicated in 
figure 1. The fields within the ftth layer and on its boundaries can be represented in matrix form as 

[S„(r)] = [a„(r)]-[C„], (55) 

where [S n (r)] is a column matrix of the tangential field components, [C n ] is a column matrix of 
coefficients or wave amplitudes and [a n (r)] is a matrix of the functions or the solution matrix. At 
the two boundaries, r n -\ and r n of the nth layer (55) may be written as 

[S n (r n -i)] = [a„(r n -i)] ■ [C n ] (56) 

[S n (r n )] = [a n (r n )]' [CJ; (57) 
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FIGURE 1. A spherically stratified medium. 



UPPER MEDIUM 




[S n (r n -i)] is related to [Sn(r n )] by eliminating [C n ] from (56) and (57) as 

[S„(r„)] = [a tl (r n )] ■ [d^rn-i)] ' [S n (r n -i)] 



or 



(58) 



[S n (r n )] = [dn] ■ [Sn(r n -i)]. (59) 

The boundary condition at r n -i between the layers (n) and {n— 1) is given as 

[S n (r n _i)] = [S n -i(r„_i)] (60) 

if there are no sources on the boundary. When combining (59) and (60), the fields at the boundary 
r n are related to the fields at the boundary r () as 

[S n (r n )] = [dn] • [dn-i] . . . [A] ■ [4] ■ [a tt (r )] • [C«] = [b n ] ■ [CJ. (61) 

The column matrix is assumed of the form 



A>) 
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(62) 



There are only upgoing waves in the upper medium and the coefficient matrix [C„] becomes 



[C H ] = 



(63) 
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The matrix [a„(r)] is obtained from (31) to (37) and (44) to (50) as 

x exp (ik+r) — x exp (— ik+r) — WG' exp (ik-r) 

y exp (ik+r) y exp (— ik+r) — zE' exp (ik-r) 

— xH' exp (ik+) xH" exp {—ik+r) —w exp (ik-r) 

yF' exp (ik+r) yF" exp (— ik+r) z exp (ik-r) 



[a n (r)] = 
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and where m= ' or ". The symbols Aj and kj 1 (with 7=+ or—) designate solutions of (23) and (24) 
respectively. After computing the inverse matrix [a~ 1 (r)], the matrix [d n ] can be computed 
explicitly. The expressions of [d n ] are relatively simple if N' =N" = N(N = E, F, G, M), as can 
be seen from appendix 1. In the limiting cases of negligible displacement currents (qj << 1) 
and of uncoupled TE and TM modes (N = 0), the matrix elements are in agreement with available 
results [Thompson, 1963] . In the ionosphere layers where TE and TM modes are coupled the 
differences between N' and TV" can be ignored and the [b n ] matrix of (61) can be computed either 
by multiplying out a sequence of [a n (r n )] and \_a~ l (r n -i)] matrices or by using the explicit ex- 
pressions for the [d n ] matrix. 
Thus 



[b n ] = [a n (r„)] • [a n l (r n -i)] ■ [b n -i] 



(73) 



[b n ] = [dn\- [bn-,]. 



(74) 



In either case small elements of [b n ] are computed by substraction of large numbers, which de- 
creases the accuracy of [b n ] after repeated multiplications in particular for frequencies in the VLF 
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range at higher altitudes where the ionospherie wavelength becomes small. This difficulty can 
be avoided by evaluating the product (74) analytically, and it is shown in appendix 2 that an element 
bij\„ of [b„] is obtained from the corresponding element of [6„ i] by a single multiplication as 



Oij\n 6ij ' bij\„- 



(75) 



provided that the differences of N(N=E, F, G, H) between the two adjacent layers are negligible. 
Starting out with [a n (r )], [b n ] is computed initially by (75). At lower heights the computa- 
tions are changed to (73), and finally one obtains the matrix [b n ] at the lower ionosphere boundary 
Yn which relates the tangential field components E$ 9 7/^, E<j>, and He at this height to the amplitudes 
of the upgoing waves A and B of the upper ionospheric layer. The elements of an impedance 
matrix which characterizes the ionospheric boundary can be deduced from (61), (62), and (63) as 



Zn = 


Ee 
He 


Zn=~ 


Eft, 


Z n = 


El 


Z-n=- 


E^ 



_ 613621 — 611623 

//</,=<) 62 1 643 — ^23643' 



n<t>=o 



. 62 1633 — 623631 
623641 — 621643' 

^ 6 11643 — 61,-1641 
//0=() 621643-6236.il 

_ 631643^633641 
He =0 62364 1 — 621 643 



(76) 

(77) 
(78) 
(79) 



where b t j designates an element of the b n matrix in (61). These impedances are related to iono- 
spheric reflection coefficients, which lead to a rather involved modal equation for determining the 
propagation characteristics of the coupled TE and TM modes in the earth-to-ionosphere waveguide 
[Wait, 1963a]. For vertical dipole excitation, TE modes are not excited and also the coupling 
between TE and TM modes has been shown to be negligible for a single layer ionosphere [Schu- 
mann, 1955; Wait, 1962, p. 269]. This approximation is assumed to apply also to a stratified 
ionosphere, and the propagation characteristics are determined from the impedance Z12 in the 
same way as for an isotropic ionosphere, although the anisotropic structure of the ionosphere is 
considered in establishing Z V1 . 

The surface impedance Z J2 is computed initially by setting S=\ in v(v+ l) = (koaS) 2 and Z v > 
is used in the appropriate modal equations of ELF or VLF waves for obtaining a first estimate of 
S = S\. Using Si, the computation oiZ\% gives S2, and the process is repeated until the sequence of 
S n converges [Galejs, 1964b]. 



4. Discussion 
4.1. Extremely Low Frequency 

For frequencies in the ELF range Roots S of the modal equation are shown in figures 2 and 3 
for the ionospheric model of electronic conductivity of figures 1 and 2 of Galejs [1964b]. In the 
isotropic ionosphere model cr remains unaltered; cri^cro and cr? approaches 0. The anisotropy 
of the ionosphere increases ReS (ReS = c/i;, c = velocity of light, v = phase velocity) and also ImS, 
that is proportional to the attenuation constant. These increases are most pronounced for the 
nighttime ionosphere model in the lower frequency range shown in figure 3. The corresponding 
attenuation constants are plotted in figure 4. 
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FIGURE 2. Roots of the modal equation, electronic day- 
time conductivity for F = 0. 



Figure 3. Roots of the modal equation, electronic night- 
time conductivity for F = 0. 




FIGURE 4. Attenuation constants, electronic conductivity 
for F=0. 



FREQUENCY - C/S 



The resonant frequencies of the earth's ionosphere cavity and its (^-factors can be determined 
from the computed values of ReS and ImS [Galejs, 1964a, 1964d]. These calculations are made 
also for ionosphere models of electronic plus ionic conductivities at/=0 and/=20 c/s (figs. 1 to 4 
of Galejs, 1964b). The calculated resonance frequencies f n are listed in table 1. The isotropic 
day or night models exhibit too high resonance frequencies, but the average for the anisotropic 
day or night models gives nearly correct results although the first resonance occurs near 7.7 c/s. 
At night the various anisotropic models give different results, while their results are nearly the 
same for daytime or for the isotropic models. The calculated (^-factors are listed in table 2. At ■ 
daytime the (^-factors are lowered slightly due to anisotropy. The losses are principally due to 
absorption in the lower D region, where the collision frequencies are high and the magnetic field 
has only small effects. At nighttime the anisotropy reduces the (^-values drastically and the ions 
have a considerable effect. For the isotropic models the energy remains below the ionosphere, 
the losses in the lower ionosphere are small, and the (lvalues are high. For the anisotropic 
models some energy penetrates the ionosphere and escapes from the earth-ionosphere cavity. 
Hence the Q figure is lowered. The present calculations assume the radial magnetic field of the 
polar regions to be extended over the whole surface of the globe and thus the anisotropy effects 
of this model should be larger than expected in nature. 
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RESONANCE FREQUENCIES OF THE SPHERICAL SHELL BETWEEN THE EARTH AND IONOSPHERE 
IN THE PRESENCE OF A RADIAL MAGNETIC FIELD 



Q-FACTORS OF EARTH- IONOSPHERE CAVITY RESONANCES 
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4.2. Very Low Frequency 

For VLF waves the calculations are made first for a homogeneous sharply bounded ionosphere 
model in order to obtain a comparison with past results that are based on the quasi-longitudinal 
approximation [Wait, 1962]. The data shown in figure 5 are calculated for the same ionosphere 
parameters as in figure 12, page 276, of Wait [1962]. The tensor components of ionospheric con- 
ductivity (3) are related to the parameter 



V 



'+(? 



-1/2 



(80) 



do = € { )CO r VI + (ix)nlv) 2 
(T^CTo/ll+icO///^) 2 ] 
Cfy = CF\Cx)ull> 



(81) 
(82) 
(83) 
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Figure 5. Attenuation rates, sharply bounded homo- 
geneous ionosphere, o> r = constant. 



FIGURE 6. Attenuation rates, sharply bounded homo- 
geneous ionosphere, o- = constant. 
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where coo is a plasma frequency, coh gyrofrequency, and v is collision frequency. There is a good 
argument with the results of Wait for frequencies above 15 kc/s, but the attenuation rates of figure 
5 are somewhat higher for frequencies near 10 kc/s. This may be due in part to the differences 
between the modal equation used by Wait [1962, eq 19, p. 275] and the equations (45) to (50) of 
Galejs [1964b] which are used for the calculations shown in figure 5. (The differences of attenua- 
tion rates that are computed using such modal equations can be seen for a different ionospheric 
model from figure 14 of Galejs [1964b].) 

The attenuation rates can be computed also for cr constant in the isotropic and anisotropic 
models, as indicated in figure 6. The attenuation rates are decreased due to the anisotropy in 
figure 6. 

Higher ratios of (u>h\v) may be appropriate for the nighttime ionosphere when VLF waves 
may reach higher ionosphere regions where collision frequencies are lower. As an example, 
calculations are made for cx)hIv= 10 and for the lower ionosphere boundary at ja = 90 km, as shown 
in figure 7. There is a significant decrease of the attenuation rates due to the magnetic field, 
particularly for lower ionospheric conductivity. 

Propagation parameters are also calculated for the profiles of electronic conductivity for 
/=0 shown in figures 1 and 2 of Galejs [1964b]. At daytime there is hardly a change due to the 
radial magnetic field, but the nighttime attenuation rates are significantly decreased. This is 
in line with the data shown in figure 7 for a sharply bounded homogeneous ionosphere model. 
There are uncertainties associated with nighttime models, and the numerical differences shown in 
figure 8 may be higher than observed experimentally. 

The attenuation rates are also calculated for a daytime ionosphere model which considers 
electronic and ionic conductivities at f= 15 kc/s [fig. 5 of Galejs, 1964b]. Selection of two dif- 
ferent boundary heights y\ of the ionosphere model (ionosphere conductivity is assumed to be 
zero for y < y N ) causes only minor changes between calculated attenuation rates, as is seen from 
figure 9. The ionospheric anisotropy decreases a for y^ = 60 km, but it has a negligible effect 
for jn = 50. The effects of anisotropy are comparable to the uncertainty of the propagation param- 
eters associated with the selection of two different boundary heights y ; y. 
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Fig URE 7. Attenuation rates, sharply bounded homo- 
geneous ionosphere, a q — constant. 



Figure 8. Attenuation rates electronic conductivity 
for F = 0. 
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Figure 9. Attenuation rates, electronic plus ionic day- 
time conductivity for F = 15 kc/s. 
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Figures 5 to 9 have shown only the attenuation rates which are proportional to ImS. The 
changes of ReS due to the radial magnetic field are negligible and the ReS figures for an isotropic 
ionosphere [Galejs, 1964c] apply also to the case of radial magnetic field. 
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5. Appendix 1. Matrix [d n ] 

After computing the inverse matrix [a~ l {r n -\)], the elements of the matrix 

[dn] = [a n (r n )] ' [a'Krn-i)] 



are computed as 



d\ i = n , Tr ^ — — [cos k+h + HG cos k-h] 
1 + HG r n 



dv>— — 
di3 = 
du = 



r n -\ 



l + EF r„ 
G r n 



X n j , , W n GF . 

sin k+n-\ sin k-h 



1+i/G r n 

i r n -i 



l+EF T n 



y„ z n 

[ — cos k+h-h cos k-h] 



x n tj . iii w n G . , L 
sin k+h-\ sin k-h 

Zn 



Yn 



j i r n -\ 

dn ~~TTHG r n 



z n EH . 



— sin k + h-\-— sin k-h 



dzz — -, , riT^ ~^~~ t cos k+h + EF cos k-h] 



d 2 * = 



l+EF r n 
i rn-j 

1 + HG V n 



sin k+h sin k-h 
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(84) 

(85) 
(86) 
(87) 
(88) 
(89) 
(90) 
(91) 



761-367 0-65-4 



diA— , F p [cos k+h — cos k-h] 

dsi = i i TTr , [ ~~ cos k+h -h cos A:_/i] 



d,v> '- 



1 + //G r M 

t r n -\ 
1 + EF r n 



sin k+h sin k-h 



dx\ — llr , [cos k-h + HG cos k+h] 

1 -r tlLr r n 

dM ~T+EF~ 



d 4 \ — 



i r n -\ 



1 + //G r n 



k+h-\ — - sin k-h 

Z n J 

— sin k+h-\--^— sin k-h 

n W n J 



X n HE . 

sin 



c?42— t i t?t? ~ — [ cos k+h — cos k-h] 



d^ — 
d 44 = 



\+EF r„ 

i r w -i 
1+#G r n 

1 r n - 
\+EF r n 



— sin k+h + — sin k-h 

%n Wn 

- [cos k-h-\-EF cos k+h] 



(92) 
(93) 
(94) 
(95) 
(96) 
(97) 

(98) 
(99) 

(100) 



where the symbols are defined as in (64) to (72) and where h = r n -i — r n . However, these expres- 
sions are valid only if N'=N" in (69) to (72). 

When neglecting ionospheric displacement currents q/s are small and 



F=H=-G=-E 



(101) 



and the matrix elements simplify to 



r„-\ 



dt \ = d>-2 — d:n = di4 = „ — (cos k+h + cos k-h) 
2r„ 

, , _ icojU-o r„-i 

"12 — «34 — q 

I r„ 



sin k+h sin k-h 



d t . 



d\A — dvz — 



k+ 



- d>>4 = d 42 = o — — [ — cos k+h + cos k-h] 



sin k+h sin k-h 



(Q/JLq r n - 
2 r n L £+ 



d2\= — d 4 3 = — 7; [&+ sin /c+/i 4- Zc_ sin k-h] 



d 2 3 = d 4 i=- 



2o)fx r n 

1 r n -\ 
2co/x r„ 



[/c+ sin k+h — k- sin A:_A]. 



(102) 
(103) 

(104) 
(105) 

(106) 

(107) 



The coefficients (102) to (107) are in agreement with results of Thompson [1963], which were 
originally obtained starting out with the field representation (38) to (43). 

For uncoupled TE and TM modes, E, H, F, and G approach zero and the matrix elements 
(85) to (100) become 
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d n = d22= JL - L cos k+h (108) 

ico/x (l — <7;0 r„_i . . , 
di2 = 7 sin k+h (109) 

rfi:j = di4 = d 2 z = diA = d:u = d.vi = d 4 i=d 4 2 = (110) 

d 2 i=- 77 r-^ sin k+h (111) 

icofioil — g 3 ) r n 

d^ = du = J1 -- cos k-h (112) 

A;- r„ 7 

A- r -i 
rf43= — ' " — sin k- h. (114) 

After expanding the sine and cosine functions and substituting the values of k+ and k- from (28) 
and (30), the coefficients (108) to (114) are in agreement with the free space results of Thompson 
[1963], which were originally obtained from a Taylor series expansion of Hankel functions repre- 
senting the fields between earth and ionosphere. The matrix elements (108) to (114) were derived 
using the planar approximation (14). These expressions remain valid in the ELF range in the 
spherical shell between the ground and ionosphere. For VLF, these expressions should be used 
only in thin shells like in the vicinity of the ionospheric boundary. 

6. Appendix 2. Multipliers e i} 

After carrying out the matrix multiplications 

[6i] = [rfi] • [a n (r )] = [di] • [b {) ] 
[4i] = [*]-[6i] 

: (lis) 

[bn] = [dn]-[b n -l] 

it is found that the ratio between an element of the matrix [b„] and the corresponding element 
of the matrix [6>,-i] is given by 

eij = bij\nlbij\ n -l (116) 

and it does not depend on other elements of the [6,*-i] matrix if E, F, £, and H of two adjacent 
layers are the same. The multipliers eij are obtained as follows: 

/ 7 , ■ X n Jn-\ • , u \ r _n-l ,_,_ 

en = e : u = I cos k+h — i sin k+h I (117) 

€vz = evz — I cos k+h + i — - — — sin k+h J — — (118) 

( .W n Zn-1 . r L \ r u -\ 

Ci3 = e:v.i = cos k-h — i sin k h (119) 

\ Z n W n -\ J r n 
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ei4 = e 3 4— ( cos k-h + i— -— ^— -sin k-h) — — (120) 

ei\ = C4i = ( cos k+h — i — — — sin k+h) J! -^ (121) 

X n Yn-1 I r n 



e 2 2 



— 642 = ( cos /c+A H- i — - — - — sin /c+A I — — (122) 

\ X n Yn-1 I r n 



e-iz — e 4 3 — ( cos k-h — i — — — sin /c_A I — — (123) 

\ W n Z n -1 ] r n 

( , . Z n Wn-1 . j L \ T n -\ 

^24 — ^44= cos k-h-\-i — sin k-h (124) 

\ W n Zn-1 I T n 
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